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SUMMARY 


The digital versions of optimal-linear-regulator theory and eigenvalue 
placement theory are applied to the Mach number control loop of the National 
Transonic Facility cryogenic wind tunnel. The control laws developed are eval- 
uated on a nonlinear simulation of the tunnel process for a typical test 
condition and are found to significantly reduce the open-loop time required to 
achieve a Mach number set point. 


INTRODUCTION 

The National Transonic Facility (NTF) , figure 1, will be a continuous- flow 
cryogenic wind tunnel designed to satisfy the basic airfoil research and develop- 
ment needs of the federal government and industry. The tunnel (refs. 1 and 2) 
will operate on the basic principle that Reynolds number is inversely propor- 
tional to temperature. It is expected to achieve an order of magnitude increase 
in Reynolds number over existing tunnels at reasonable values of dynamic pres- 
sure by injecting liquid nitrogen as a coolant to maintain the test gas at cryo- 
genic temperatures. The NTF will be controlled by a digital canputer in order 
to achieve high productivity and to implement control of dynamic pressure, Mach 
number, and Reynolds number, or other equivalent sets of process variables. 

A plan view of the tunnel control circuit is shown in figure 2. The NTF will 
be a closed-circuit fan-driven pressure tunnel capable of total pressures up 
to 896 31 8 Pa (130 psi). It will have a 6.25-m^ (67.24-ft^) slotted test section 
which will be 61 m (200 ft) long and 14.6 m (48 ft) wide. It will operate from 
room temperature to 88.70 K (-300^F) . Liquid nitrogen (LN 2 ) will be sprayed into 
the circuit upstream of the fan (a) for the initial cooldown phase, (b) to 
balance the heat rise associated with the gas compression by the fan, (c) during 
steady-state temperature regulation, and (d) to act as the primary mode of tem- 
perature control in the test section. During the ambient temperature operation 
(air mode) of the NTF cooling functions, (b) , (c), and (d) are taken over by a 
heat exchanger, located in the stagnation section, and an external cooling tower. 
No LN 2 is consumed in the air mode. Pressure control is provided by venting gas- 
eous nitrogen (GN 2 ) from the tunnel to the outdoors. Coarse Mach number control 
will be provided by fan motor speed regulation, and fine vernier Mach number 
variations will be controlled with small changes about a constant reference guide 
vane angle. 

Digital computers and microprocessors will be employed for autcxnation and 
process control of the NTF. Among the functions planned for computer operation 
are supervisory set-point programming, monitoring process limits, sequencing and 
timing, and automation of the start-up and shutdown operational sequence. Addi- 
tionally, the digital capability will offer the opportunity to apply multivari- 
able control techniques to the wind-tunnel process. The purpose of this report 
is to demonstrate the use of two multivariable methods to design control laws for 
a particular NTF operating condition. 



Analytical studies of the NTF process have produced a transfer function 
describing, in the s-domain, the dominant interaction of the test-section Mach 
number and fan guide vane angle when the process is perturbed from a steady- 
state operating condition. This transfer function is employed herein to 
demonstrate the application of multivariable design methods to the Mach number 
control loop of the NTF. Initially, in order to overcome a design difficulty 
presented by a transport lag appearing in the transfer function and to gen- 
erate time-domain equations needed for multivariable design, finite-difference 
equations are derived which approximate the Mach number to guide vane angle 
characteristics for zero-order hold guide vane actuator and actuator rate 
inputs. Optimal linear regulator and eigenvalue placement techniques (refs. 3 
and 4) are next applied to derive multivariable feedback control laws for Mach 
number regulation about a typical steady-state operating condition. System 
responses are calculated from the linear finite-difference design equations and 
a nonlinear simulation, which solves the Navier-Stokes partial differential 
equations for unsteady flow to represent the fluid flow dynamics. 


SYMBOLS and ABBREVIATIONS 

A response matrix for state equation in variable x (see eqs. (A4) , 

(31) , and (34)) 

Aq,Ai defined by equations (15) and (16) 

yN ^ /N 

Aq,A-|,A 2 defined by equations (20), (21), and (22) 

B control effectiveness matrix for state equation in variable x (see 

eqs. (A4), (32), and (35)) 

Bq defined by equation (17) 

B-| defined by equation (A1 ) 

Bq defined by equation (23) 

ys 

B] defined by equation (A13) 

col column matrix 

diag diagonal matrix 

f feedback gain 

GN 2 gaseous nitrogen 

h time delay, sec 

i nonnegative integer 

j =Npr 
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k gain in transfer function defined by equation (1), deg“^ 

LN 2 liquid nitrogen 

M Mach number 

Mf set-point Mach number 

m natural number defined by equation (18) 

N positive integer, such that NA = T 

P defined by equation (45) 

Pt,f set-point total pressure, Pa 

Q weighting matrix in equation (45) 

rpm revolutions per minute 

s Laplace transform parameter 

T final time for matrix equation (A1 ) , sec 

Tt,f set-point total temperature, K 

t time, sec 

tg settling time, sec 

u control variable 

X state vector defined by equations (33) and (36) 

X state vector defined by equation (14) 

2 defined by equation (24) 

A sampling period, sec 

6M Mach number perturbation 

60 guide vane angle perturbation 

60^ guide vane angle actuator-input perturbation 

C damping ratio 

0 guide vane angle, deg 

0;^ guide vane angle actuator input 
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6f set-point guide vane angle, deg 

T time constant, sec 

(}) matrix defined by equation (46) 

(jL) frequency, rad/sec 

tera,b) a^t<b 
« approx imately equal 

A dot over a symbol denotes differentiation with respect to time. 

A prime after a symbol denotes a matrix transpose. 

MACH NUMBER TRANSFER FUNCTION CHARACTERISTICS 

This section describes some basic characteristics of the NTF Mach number 
control loop. The information to be presented was obtained by means of private 
communication from George Gumas, Pennsylvania State University, Middletown, 
Pennsylvania. 

Consider the NTF process to be operating in a steady-state condition such 
that, for some fixed fan speed, constant values of guide vane angle, liquid- 
nitrogen injection rate, and gaseous-nitrogen vent rate are chosen to maintain 
all tunnel states at zero time rate of change. When the guide vane angle is 
perturbed from its steady-state value, while maintaining the liquid-nitrogen 
injection and gaseous-nitrogen vent rates, the transfer function between test- 
section Mach number and fan guide vane angle closely approximates the following 
equation if no acoustic modes are excited and perturbations remain small: 


6M(s) ke 

60 (s) (Ts + 1) 


( 1 ) 


The symbols 6 M(s) and 60 (s) represent, in the s-domain, the perturbations 
away from steady-state Mach number Mf and guide vane angle 6f and are given 
as 


6M(s) = M(s) - Mf/s 


( 2 ) 


and 


50 (s) = 6 (s) - 0f/s 


(3) 
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The constants k and T in equation (1) depend parametrically on the operating 
condition. Tables I to III (provided by George Gumas) show k and T values 
for steady-state test-section total temperature and Mach number Mf. The 

transport lag h is approximated by 


h = 4.29/\pP^ 


(4) 


Additionally, the guide vane angle is set through an actuator with the 
transfer function approximated by 


60 (s) 
60^(s) 


= w2/(s2 + 2 ^ 0 ) s + 0 ) 2 ) 


(5) 


where 


0) =6.0 rad/sec 


and 


C = 0.8 


and where 


59a(s) = 9 a(s) - 9 f/s 


( 6 ) 


and 0;^ is the transient value of the input to the guide vane actuator. 

Operational constraints limit the values of guide vane angle 6 and guide 
vane angle rate 6 to lie within ±30® and ±30®/sec, respectively. 

Control laws are derived for Sd and 56;^ which regulate the Mach number 
about the steady-state value Mf. Block diagrams for the composite transfer- 
function sequences are shown in figure 3. The multivariable design methods to 
be applied require constant-coefficient time-domain equations for equations (1 ) 
and (5) which are developed in the section which follows. 
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MATHEMATICAL MODELS FOR MULTIVARIABLE CONTROL DESIGN 


In this section, time-domain vector-matrix equations of the form 


xF (i + 1)A] = A x(iA) + B u(iA) 


(i = 0,1 ,...) (7) 


are developed from the dynamics of the test-section Mach number expressed in the 
s-domain by equations (1 ) and (5) . The matrices A and B in equation (7) are 
constant and A is the sampling period for applying either of the zero-order 
hold controls; 

guide vane actuator input 


u(iA) = 56j|^{iA) 


(8) 


or guide vane actuator rate input 


u(iA) = 60^ (iA) 


(9) 


The Mach number response 6M(iA) is one of the elements of the state vector x. 
The delay equation 


T 6M(t) + 6M(t) = k 69{t - h) 


( 10 ) 


and the ordinary differential equation 

60 (t) + 2C0) 50 (t) + 0)2 50 (t) = (Jj2 60^ (t) 


( 11 ) 


have the same input-output characteristics as equations (1) and (5), respec- 
tively. Equations (10) and (11) can be written in state-variable form as 

x{t) = Aq x(t) + A-| x(t - h) + Bq u(t) (12) 


where 


u(t) = 60A(t) 


(13) 
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X(t) 


( 14 ) 


Ao 



Xl (t) 


6M(t) 

= 

X2(t) 

= 

66 (t) 


xa(t) 


60 (t) 

f-l/'T 

0 



0 0 1 
2 

0 -W -2 


0 


A) = 


0 


0 


k/T 0 

0 0 

0 0 


and 


0 


Bo = 


0 


032 


( 15 ) 


( 16 ) 


( 17 ) 


Multivariable techniques exist which could be applied directly to equa- 
tion (12), but the calculations are greatly complicated by the presence of the 
time delay h. One approach to overcoming the problem of time delays is to 
use a zero-order hold structure for the control and approximate the time delay 
by an integral multiple of the sampling period; that is. 


nA » h 


( 18 ) 


for some nonnegative integer m. 

This approach to treating time delays is described in the appendix and, when 
applied to equation (12), produces 

x[ (i + 1 )A] = Aq x(iA) + A-| x[ (i - m)A] + A 2 xf (i - m + 1 )A] + Bq u(iA) (1 9) 
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where 


e'VT 


In = - I 


coo \jl - y/>Ji" _' ^2je-CAo» ^Aw \jl - j 

(“to y/\|l - qZje-cAiiJ sin ^to\|l - ^2j 


^l/w\|l - ^je-tAo) sin ^Aaj\|l - 
-CAd) cos ^Au)\jl - ;2j _ ^5y/\|l - 52je“CAu> sin ^AuNJi - 


( 20 ) 


Al 



0 

0 


2T 

0 


0 

0 


0 0 


0 


( 21 ) 


A2 



0 

0 


0 


Ak 

2T 

0 

0 


0 

0 

0 


and 






eAoT 


dT Bq 


(22) 


(23) 


Because of the sparseness of A-| and A2, the delayed terms in equation (19) 
involve only the variable X2« After introducing m (eq. (18)) additional vari- 
ables and equations of the form 


z-|[(i + 1)A] = X 2 f(i - m + 1 ) A] = Z 2 (iA) 
Z 2 f (i + 1 ) A] = X 2 [ (i ~ m + 2 ) A] = 23 (iA) 

J ) 

Zm_i r ( i + 1 ) A] = X 2 r ( i - DA] = Zm ( 
z„,r (i + 1 ) A] = X 2 (iA) 


( 24 ) 
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then equations (19) to (24) can be combined to construct an equation of the form 
of equation (7) with state vector 


x(iA) 


xi.(iA) 


6M(iA) 

X2(iA) 


60 (iA) 

X3(iA) 


60 (iA) 

zi (iA) 


60 r (i - m)A] 

Z2(iA) 

m 


60 r (i - m + 1 )A] 

• 

• 

• 

Zm(iA) 


• 

• 

60f (i - 1)A] 


(25) 


and guide vane actuator input control given by equation (8) . 

The coefficient matrices depend parametrically on the operating condition 
and sampling period chosen. A typical set of conditions for NTF operation within 
a linear region of the data sets of tables I to III was selected for this study. 
These conditions are as follows: 


Mf = 0.9 

Tt,f = 166.67 K (-160°F) ) 

Pt,f = 517 107 Pa (75 psi) 

’ J 

The fan is operated at synchronous speed. Table II lists 


(26) 


T = 1 .964 sec 
k = -0.0117 deg-'' 


(27) 


and equation (4) yields 


h = 0.33 sec 


(28) 


9 



Then, for a sampling period of 


A = 0.1 sec 


(29) 


the selection of 


m = 4 (30) 


causes the actual value of iA + h to fall within the fiA - mA, iA - (m - 1)A] 
interval in equation (19). For the conditions shown in equations (26), the 
set-point guide vane angle 0f is 1 .93°. Using the foregoing numerical data 
and the guide vane actuator as the control u(iA), it can be determined that 


0.95035 

0 

0 

-0.002831 

-0.0002779 

0 

0 

0 

0.86976 

0.06055 

0 

0 

0 

0 

0 

-2.1798 

0.28848 

0 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

0 

1 

0 

1 

0 

0 

0 

0 

0 


(31) 


and 


with 


B = 


0 

0.1302 

2.1798 

0 

0 

0 

0 


(32) 
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0 


( 33 ) 


the control u(iA), a similar proce- 


- 0.0002831 

0 

0 

0 

0 

0 

0 

0 


- 0.0002979 

0 

0 

0 

1 

0 

0 

0 



( 35 ) 


n 



with 


x(iA) 


6M(iA) 

60 (iA) 

60 (iA) 

60^ (iA) 

60 [ (i - 4) A] 
60[(i - 3)A] 
60 [(i - 2) A] 
60 [(i - 1)A] 


(36) 


Now that equations (1) and (5) have been combined to form equation (7), the 
multivariable design techniques are applied in the section which follows. 


DESIGN AND LINEAR SIMULATION 

Equation (7) , with coefficients and state variables defined either by equa- 
tions (31), (32), and (33), or (34), (35), and (36), represents variational equa- 
tions linearized about the set point of equation (26) . The design problem 
considered in this report is to assume that the system is initially perturbed 
from the operating point and derive, for either equation (8) or equation (9), 
constant-gain state-variable-feedback control laws which cause the Mach number 
and guide vane values to return rapidly to the operational set without violating 
the constraints 


-30° ^ 0(t) ^ 30° 

and 

-30°/sec i 0(t) ^ 30°/sec 
In other words, control laws of the form 
u(iA) = -f x(iA) 
are sought such that 


(37) 


(38) 


(39) 
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lim x(iA) •* 0 

2j>co 


( 40 ) 


without violating contraints (37) and (38) when 


x[(i + 1)A] = [a - Bf] x(iA) 


(41) 


and 


x(0) fi 0 


(42) 


Two techniques, the optimal-linear-regulator and eigenvalue placement algorithms 
(refs. 3 and 4) are employed to form the control gains. 

The optimal-linear-regulator algorithm finds the feedback gain by selecting 
the controller which minimizes a quadratic performance index 


J 


E 

n =n 


{x' f (i + 1 ) A] Qx[(i + 1)A] + u^(iA)} 


(43) 


subject to equation (7.) , where Q is a nonnegative-definite symmetric weighting 
matrix. When certain controllability conditions (ref. 4) are satisfied, the gain 
is given by 

f' = (1 + B’PB)"‘'b'PA (44) 


with 


p = (|>'p<t) + ff' + Q 


(45) 


and 


(j) = A - Bf 


(46) 


Software for solving equations (44), (45), and (46) is described in reference 3. 
Parametric studies of the elements of Q provided the choice of 


Q = diag (10^,0, ... ,0) 
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as a weighting matrix yielding control gains which, upon simulation, typically 
produced responses satisfying the constraints -30° S 0 (t) S 30° and 
-30°/sec ^ 6 (t) S 30°/sec. For A, B, and x as given in equations (31), 
(32) , and (33) , the optimal regulator procedure gave the control law 

66;^(iA) = 31 .73 6M(iA) - 0.07685 66 (iA) 

- 0.006603 60 (iA) - 0.009451 66 f (i - 4)A] 

- 0.01989 60f (i - 3)A] - 0.02093 66T (i - 2)A] 

- 0.02202 66f (i - 1)A] (48) 


For A, B, and x as given in equations (34) , (35) , and (36) , the optimal 
regulator procedure gave the control law 


66\(iA) = 27.31 6M(iA) - 0.06776 66 (iA) 

- 0.006142 66* (iA) - 0.6808 60^^(iA) 

- 0.008148 60[ (i - 4)A] - 0.01712 66 f (i - 3)A] 

- 0.01801 66f (i - 2)A] - 0.01895 50 f (i - 1 )A] (49) 


In the eigenvalue placement algorithm (refs. 3 and 4 ) , the designer selects 
the eigenvalues of (A - Bf*) which in turn determines the speed of closed-loop 
response. It was generally found that significant improvement in the Mach number 
response speed could not be obtained using control law (8) and the eigenvalue 
placement algorithm without violating constraints (37) and (38). Results are 
therefore only presented for eigenvalue placement using the actuator rate input 
(eq. (9)). The eigenvalues of A given by equation (34) are as follows: 


Eigenvalue 
0.950 
1 .0 

0.579 ± j(0.218) 
0 
0 
0 
0 


Mode 

Mach number 

Integrator 

Actuator 

> Delayed states 

J 
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The delayed states are in deadbeat form and need no alteration. The actuator 
dynamics cannot be physically altered. Only the Mach number and integrator modes 
are available for design. The best results were obtained when f was selected 
to cause both these eigenvalues to take on the value 0.9. The resulting control 
law is 


66j^(iA) = 31 .833 6M(iA) - 0.0789 60 (iA) 

- 0.007218 60(iA) - 1.490 50^(iA) 

- 0.009482 60[(i - 4)A] - 0.01995 60f(i - 3)A] 

- 0.02100 6er(i - 2)A] - 0.02209 60[(i - 1)A] (50) 


In order to demonstrate the characteristics of the optimal regulator and 
eigenvalue placement control laws, the transient response from a steady-state 
condition at Mach 0.8 is calculated. Initially, 


6M(0) = M(0) - Mf 
= 0.8 - 0.9 

= -0.1 (51) 


Moreover, the steady-state condition implies that, for t S O, 

"X 

6M(t) = -0.1 
60 (t) =0 \ 

60^(t) = 60 (t) 

J 

Then, from equation (31) or equation (34), it can be seen that 


6M[(i + 1)A] = 0.95035 6M(iA) - 0.0002831 60f(i - 4)A] 
- 0,0002979 60[(i - 3)A] 


(52) 


( 53 ) 


Substituting 

60[(i - 4)A] = 60[(i - 3)A] = 60(t) 
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and 


5Mf(i + 1 ) A] = 6M(iA) = -0.1 


into equation (53) gives 


60(t) = 8.55° 


( 54 ) 


for t ^ 0. 

The open-loop Mach number and guide vane angle responses to initial condi- 
tions of equations (51) and (52) using control inputs 


6^(t) = ^(t) - 6f = 0 


(55) 


are shown in figure 4. There it is shown that the variational Mach number dynam- 
ics^ linearized about Mf = 0.9^ is self-regulating; that is, 6M(t) converges 
monotonically without overshoot to the origin (M = 0.9) with no 60;^ input. 

Since the eigenvalues of equation (31) are within the unit circle of the complex 
plane, self-regulation occurs for all initial conditions 6M(0) . Convergence to 
0.9 is assumed to have occurred when a settling time tg is reached such that 

!6M(t)j ^ 0.002 for t ^ tg (56) 

For figure 4, settling time tg is 8.4 sec. 

Figure 5 shows the closed-loop Mach number and guide vane angle responses 
with the optimal regulator actuator position input control law (48) . The Mach 
number response monotonically increases to Mf with a settling time of 5.7 sec 
or an improvement of 2.7 sec over the open- loop results. Guide vane action is 
smooth with a steady-state value of 0£ = 1 .93^. The guide vane angle initially 
decreases in order to increase Mach number. 

Figure 6 shows the responses with the optimal regulator actuator rate input 
control law (49) . Mach number convergence is no longer mono tonic. Overshoot 
(M > 0.9) occurs at about 5.4 sec and does not settle down to within ±0.002 until 
7.4 sec. Guide vane angle convergence to Qf is slower than that shown in 
figure 4. 

The responses with the eigenvalue placement actuator rate input control 
law (50) are shown in figure 7. Mach number convergence is monotonic with no 
overshoot. The settling time is about 6.3 sec. 
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Based on the linear simulation/ control law (48) has the better 
performance - monotonic behavior with lower settling time. Control law (49) 
controls only slightly better than the open-loop response. The overshoot 
causes a time-consuming delay in convergence. Control law (50) has excellent 
characteristics but a longer settling time than control law (48). However, 
control law (50) is found to perform better than control law (48) when tested 
on the more realistic nonlinear simulation described in the next section. 


NONLINEAR SIMULATION 

The physics of actual fluid flow in the NTF produce interactions between the 
temperature, pressure, and Mach number control loops which were not accounted 
for in the previous linear simulation. In order to test control laws (48) 
to (50) in a more realistic situation, they are applied to a one-dimensional dis- 
tributed parameter model of the NTF which is based on the Navier-Stokes equations 
for unsteady flow. Control laws developed by Gumas for temperature and pressure 
regulation are employed in this nonlinear simulation to maintain total tem- 
perature and total pressure at steady flow values during Mach number transitions 
effected by control laws (48) to (50). The performance of these control laws 
when applied to the distributed parameter model is compared to that observed 
using the linear model. 

The distributed parameter model of the NTF is realized by solving the 
Navier-Stokes partial differential equations of fluid flow (including viscous 
and thermal terms) by means of McCormack's predictor-corrector method and the 
use of an artificial viscosity technique for treating shocks. Tunnel geometry, 
liquid nitrogen injection, gas venting, and fan dynamics are introduced in the 
boundary conditions. Actuator dynamics, valve calibrations, fan-flow and thermal 
relationships, heat transfer of gas to tunnel liner, flow dynamics of test sec- 
tion to plenum, tunnel-wall flow losses, screen and turning-vane losses, strut- 
sting drag losses, and sensor dynamics are all included in the model. The model 
is simulated on a vector-processing digital computer requiring approximately 
11 sec of machine time for each second of simulated time. 

Gumas* control laws for total pressure and temperature act independently 
of each other and the Mach control loop. Both control loops are of the 
proportional- integral (PI) form designed for both steady-state regulation and 
state transition control. 

In the pressure control loop, total pressure is sensed in the stagnation 
section of the NTF and subtracted from the pressure set-point value to produce 
an input error signal to a PI controller. The PI controller output is summed 
with a feed-forward signal, which is a filtered value of the estimated vent- 
valve position necessary to achieve the commanded pressure set point. This 
forms the position command to the vent-valve actuator. Vent flow then estab- 
lishes total pressure. 

The temperature control loop is somewhat more complex. In the regulation 
mode, total temperature, also sensed in the NTF stagnation section, is sub- 
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tracted from a fixed- temperature set-point value, and the error is input to a 
PI controller. The PI controller output is combined with fan operating param- 
eters in a predictive control law which estimates the LN 2 flow rate necessary 
to maintain an energy balance in the tunnel. The estimated LN 2 flow rate, after 
passing through a lead-lag network which compensates for temperature sensor time 
lags, is input to an LN 2 valve control loop. This loop establishes the valve 
position needed to produce the required LN 2 flow rate. This maintains the total 
temperature. During transition control, the PI controller is inactive to prevent 
reset windup, and a slowly varying temperature set point, programmed to assure 
a smooth temperature change, is furnished to the predictive control law. The 
transition control mode was not employed in the tests reported herein. 

Figure 8 shows the open-loop Mach number and guide vane angle responses of 
the distributed parameter model, which corresponds to the linear case shown 
in figure 4. Note that the initial value for 9(0) in figure 8 differs from 
the value in figure 4. This is because the value in figure 8 is the actual 
steady-state value at Mach 0.8 for the distributed model, whereas the value for 
figure 4 is the predicted value from the linear equations. Figures 9, 10, 
and 11 show Mach number and guide vane angle responses of the distributed param- 
eter model for optimal regulator control laws (48), (49) and the eigenvalue 
placement control law (50) , respectively. The results are seen to be similar to 
the linear cases shown in figures 5 to 7- Control law (50) gives the smoothest 
and most rapid performance of the three cases shown. A comparison of settling 
times for the linear and nonlinear cases for open-loop and control laws (48) 
to (50) is presented in table IV. 


CONCLUDING REMARKS 

A transfer function for the Mach number to guide vane angle control loop has 
been employed to demonstrate the application of multivariable design methods to 
the control of the National Transonic Facility (NTF) . For a typical NTF operat- 
ing point, digital forms of the optimal linear regulator and eigenvalue techni- 
ques were applied to produce zero-order hold control laws for the guide-vane 
actuator and actuator rate inputs. When evaluated in a realistic simulation of 
the NTF process, the optimal regulator actuator input control laws and eigen- 
value placement actuator rate input control laws reduced the settling time for 
an open-loop input by more than 27 percent. However, the control law derived 
from eigenvalue placement theory gave a preferable Mach number response. Given 
a transfer matrix (relating, for example, total temperature, total pressure, 

Mach number to injection rate, vent rate, guide vane angle), the procedure 
described in this paper could be applied to produce digital control laws which 
would regulate the state variables about a given design point and account for 
cross-coupling between the modes. 


Langley Research Center 

National Aeronautics and Space Administration 
Hampton, VA 23665 
June 23, 1981 
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APPENDIX 


APPROXIMATION OF TIME-DELAY ORDINARY DIFFERENTIAL SYSTEMS BY 
FINITE-DIFFERENCE EQUATIONS 

Given a time invariant system with state x(t), control u(t), and a delay 
h > 0 such that 


x(t) = Xo(t) , u(t) = Uo(t) for -h ^ t ^ 0 


where the subscript o denotes initial condition and 


d 

— x(t) = Aq x(t) + Ai x(t - h) + Bq u(t) + B-| u(t - h) (AT ) 


for 0 ^ t ^ T, a method is presented for approximating equation (Al ) by the 
finite-difference equation 


X [ ( i + 1 ) A] = Aq X { i A) + A-| X r ( i - m) A] + A 2 x [ ( i - m + 1 ) A] 

+ Bq u(iA) + Bt u[ ( i - m)A] (i = 0,1,..., N; nA = t) (A2) 


In equation (A2) ^ Aq , , h 2 t Bq, and B'| are constant matrices^ the 

control u(t) has been assumed to have zero-order hold structure with sampling 
period Is, and m is a nonnegative integer chosen such that mA approximates 
the time delay h. Using equation (A2) and a new state vector generally struc- 
tured as 


x(iA) = col{x(iA)^ x[(i - 1 ) A] , xf (i - 2) A] , . . . ,x [ (i - m)A] , 

u[(i - 1)A], u[ (i - 2)A] ,...,u[(i - m)A]} (A3) 


the system dynamics can be stated in the form 

x[(i + 1)A] = A x(iA) + B u(iA) (i = 0,1,... ,N) (A4) 

to which modern control theory can be readily applied. 
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APPENDIX 


Begin by partitioning the interval [0 ,t] into segments tj^, (i = 0,1/ 
such that 


0 = to < t-j < t2 < ... < tu-i < tu = T 


The solution of equation (Al ) at points t^ is 


•'ti 


- h) dr 


-1 

e*0"= 

•'f • 


h■^ u(T - h) dT 


or, if A ■" 


x(ti+1 ) = ^ 


rA 

■^O^i x(ti) + Al x(ti+i - T - h) dr 


C^i 

+ 1 e^O^ Bo u(ti+i - T) dr 

*r\ 




e^O'’^ Bi u(ti+l - - h) dx 


Let ti = iA with A constant and assume 


h ®> mA 


.../N) 

(A5) 


(A6) 


(A7) 

(A8) 
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whereby equation (A7) becomes 


APPENDIX 


xf(i + 1) A] 


gAgA ^ 

•'o 

* r e*0^ 
•'0 


Ai x[(A - T) + (i - m)A] dT 


Bq u[(i + 1) A - t] dT 


Bt u[ (A - t ) + (i - m)A] dr 



(A9) 


Next impose a zero-order hold form for the control function u(t) such that 

u(t) = u(iA) = Constant, tefiA, (i + i)A) (i = 0,1,. ..,N - 1) (AlO) 
Substituting equation (A10) into equation (A9) gives 


xf (i + 1)A] = x(iA) + J' e^O^ A-j x[ (A - x ) + (i - m)A] dr 




e^O^ Bq dx 1 u(iA) + 


r 

0 


e^O^ Bi dx 


u[ (i - ra)Al (Al 1 ) 


Finally, applying the trapezoidal integration rule to approximate the first 
integral in equation (Al 1 ) 


y '*A 

gAoT 

0 


A-| xf (A - X) + (i 

- m) A] drc 



fV 

^ AqA 
- e ^ Al 
2 

x[ (i - m)A] + 

A 

i''- 


x[ (i + 1 - m)A] 


(A12) 
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APPENDIX 


equation (All) appears in the form of equation (A2) with 





/o'' 




Al 



A2 





Bo 






eAo^ 


e^O'" 


dT Bo 


dT 

J 


(A13) 
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TABLE I.- VALUES OF TIME CONSTANT T AND GAIN k 


FOR T4- = m.n K (-260OF) 


k, deg ^ 

-O, 00344 
-.00435 
-.00734 
-.00834 
-.00978 
-.01 41 5 
-.01 44 
-.01 02 
-.00836 


1 .17 
1 .1 5 
1 .69 
2.08 
2.39 
2.58 
2.76 
1 .97 
1 .69 



TABLE II.- VALUES OF TIME CONSTANT T AND GAIN k 


FOR Tt,f = 1 66.67 K (-1 600F) 







TABLE III.- VALUES OP TIME CONSTANT T 


AND GAIN k 


FOR Tt,f = 222.22 K (-60° F) 


Mf 

rpm 

T, sec 

k, deg~^ 

0.4 

240 

0.739 

-0.00290 

.5 

240 

.830 

-.00403 

.6 

360 

1 .1 50 

-.00820 

.7 

360 

1 .470 

-.00920 

.8 

360 

1 .640 

-.01 01 0 

.9 

360 

1 .690 

-.0111 0 

1 .0 

360 

1 .760 

-.01 330 


TABLE IV.- MACH NUMBER SETTLING TIMES FOR LINEAR 
AND NCa«LINEAR SIMULATIONS 


Control law 

Settling 

IN 

times, sec 
$ 0.002 

Linear 

Nonlinear 

e^ct) = 6f 

8.4 

5.5 

Equation (48) 

5.7 

4.0 

Equation (49) 

7.4 

7.0 

Equation (50) 

6.3 

4.0 







Figure 3.- Block diagrams for NTF Mach number control loop 
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jure 4.- Open-loop response - linear simulation. 
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Figure 5.- Optimal regulator for actuator rate input control law (48) 

linear simulation. 
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Jtimal regulator for actuator rate input control law (49) - 
linear simulation. 
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Figure 7.- Eigenvalue placement for actuator rate input control law (50) 

1 i near s imul a t i on . 
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